LDA+Gutzwiller Method for Correlated Electron Systems 
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Combining the density functional theory (DFT) and the Gutzwiller variational approach, a 
LDA+Gutzwiller method is developed to treat the correlated electron systems from ab-initio. All 
variational parameters are self-consistently determined from total energy minimization. The method 
is computationally cheaper, yet the quasi-particle spectrum is well described through kinetic energy 
renormalization. It can be applied equally to the systems from weakly correlated metals to strongly 
correlated insulators. The calculated results for SrV03, Fe, Ni and NiO, show dramatic improvement 
over LDA and LDA+U. 
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Despite of the successful stories of local-density- 
approximation (LDA) for simple metals and band insula- 
tors, the applications of first-principles calculations based 
on density functional theory (DFT) are limited for cor- 
related systems, due to insufficient treatment of electron 
correlation. It has long been a challenge to have a sim- 
ple scheme which keeps the accuracy and parameter-free 
character of DFT-type calculations, but can be equally 
applied to correlated and non-correlated systems. Along 
this line, significant efforts have been made in recent 
years to take into account explicitly the strong e-e in- 
teraction [H,[i [3]- Among them, the LDA+U method fl 
adopted the orbital dependent Hartree-like scheme, and 
has been successfully applied to predict the insulating 
feature of strongly correlated systems with long-range 
ordering. However, it fails for intermediately correlated 
metallic systems. The recently developed LDA+DMFT 
method 3] follow the same essence of LDA+U scheme 
(i.e, the same effective Hamiltonian), however, the inter- 
action is treated by the proper calculation of frequency- 
dependent electron self-energy through the dynamical 
mean field theory (DMFT). The method is quite success- 
ful [H, but the frequency-dependency make the method 
very expensive in practice. It is therefore important to 
develop a new method which covers from the correlated 
metal to the Mott insulator and at the same time is com- 
putationally as cheap as LDA+U. 

We will present in this paper that the Gutzwiller ap- 
proach originally proposed to study the ferromag- 
netism in correlated band and later widely used as analyt- 
ical approach for correlated model Hamiltonian Q, can 
be naturally combined with the DFT formalism based 
on the variational principal. The idea is strongly mo- 
tivated by the recent progresses for the Gutzwiller ap- 
proach: (1) the original single band Gutzwiller method 
has been generalized to multi-band case, and is exact 
for infinite lattice dimensions 0; (2) the Gutzwiller ap- 
proach can naturally cover both the noninteracting and 
the atomic limit Q; (3) the post-LDA Gutzwiller treat- 



ment (starting from the effective Hamiltonian extracted 
from LDA band) has been successfully applied to some 
systems d, 0, [l(J. Having those knowledge, we will 
show in the present paper that a fully self-consistent 
LDA+Gutzwiller method can be developed, and it is easy 
to be implemented in existing code and computationally 
simple. Implementing the present scheme in our pseudo- 
potential plane-wave code, we have calculated the non- 
magnetic metal SrV03, magnetic metal Fe, Ni, and an- 
tiferromagnetic insulator NiO. The results demonstrate 
significant improvement over LDA and LDA+U. 

The basic idea of DFT and the Kohn-Sham equation 
is to find a reference system (usually non-interacting sys- 
tem with single slater determinant wave function l^o))- 
which has the same charge density but the kinetic energy 
can be explicitly computed. In such a way, all unknown 
(both kinetic and potential) parts are moved into the 
exchange-correlation energy E xc . Unfortunately, the ex- 
act functional form of E xc is not known, and the LDA 
for E xc is not sufficient for the materials containing the 
partly filled narrow bands, such as the 3d or 4/ states. 
To improve this, the idea of LDA+U or LDA+DMFT 
method is to correct the LDA Hamiltonian H LDA by in- 
cluding explicitly e-e interaction and remove the double 
counting part from the LDA. Therefore in both LDA+U 
and LDA+DMFT, the Hamiltonian reads, 



H 



LDA+U 



H 



LDA 



I' 

ia{3 



(1) 



where \ia) are a set of local spin-orbitals with occupation 
number rii a for lattice site i, and U a p gives the inter- 
action strength between the orbitals, Hdc is the double 
couting term. 

Now the problem to be solved is to minimize the total 
energy through variational principal. In LDA+U, this is 
done with respect to the single Slatter determinent trial 
wave function l^oK which can only improve the interac- 
tion part of the total energy and leave the kinetic part 
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unchanged. In LDA+DMFT, however, both the single 
particle wave functions and the spectral density of the 
local orbitals are considered to be the variational param- 
eters. The enlargement of the variational space makes 
the LDA+DMFT work quite well for both interaction 
part and kinetic part of the total energy. But the price 
of LDA+DMFT is that we have to deal with the fre- 
quency dependence of the self energy for each local or- 
bitals, which makes it very expensive. 

An alternative way for the problem is to improve the 
trial wave-function such that the kinetic part can be 
treated as better as possible. This is what we will done 
in the following. It is known that |^o) is n °t good, 
and a better candidate for the trial wave-function is the 
Gutzwiller wave-function which is basically the lat- 
tice version of the Jastrow type wave function for quan- 
tum liquid and has been extensively used in the strongly 
correlated systems @. It is defined as 
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The projection operator P is used to reduce the weight 
of configurations for the local orbitals with relatively high 
energy, and it reads 



P = 



n[i + (A ir -i)|i,r) (i,r[ 



(3) 



where T denotes the T-th many-body configuration for 
the local orbitals in given unit cell i and A^r are the vari- 
ational parameters describing the weight for given con- 
figuration r. Therefore our task is to minimize the total 
energy with respect to the Gutzwiller wave-function: 
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Thanks to the recent progresses for Gutzwiller ap- 
proach, although the exact solution of the above expecta- 
tion is unknown, the Gutzwiller approximation (GA) is a 
very tractable approximation for the multiband Hubbard 
model and has been proven to be exact for the infinite 
dimension. In the present paper, we apply the GA di- 
rectly to the LDA Hamiltonian rather than tight binding 
form extracted from it. To do so, we only need a set of 
local orbitals \ia) where we implement the local corre- 
lation effect. Following reference Q, we can derive the 
total energy under GA as, 
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with the Gutzwiller effective Hamtonian in momentum 
(k) space given as 



H L ef DA + G = [Q 



a k 



(1 - R)]H LDA [Q + (1 - R)] 
(1 - ql)e° a \ak) {ak\ (6) 



Q = q a \ak) (ak\ , R = \ak) (ak\ 



(7) 



ak 



a k 



where Er is the eigen value and m r is the weight of the 
T-th configuration, e° is the on-site energy of each lo- 
cal orbital. The key quantities in our formalism are q a 
(0 < q a < 1), which describe the kinetic energy renormal- 
ization of local orbitals, and can be expressed in terms of 
mr through 
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The main strategy of above formalism is that, kinetic 
energy of correlated orbitals are renormalized by factor 
q Q , while all other non-correlated orbitals are taken into 
account by (1 — R) in the complete local basis projection. 
Therefore under GA, the total energy can be expressed 
again as the functional of the non-interacting wave func- 
tion but with additional variational parameters of the 
weighting factors mr, along with the following necessary 
constraint 

E < r l C « C « l r > ™r = na = E (*o| Cl k C ak |*o) (9) 
r k 

the original Gutzwiller variational parameters A^r can be 
obtained by A r = with m r = (* | |r). 

Having above equations, it is easy to obtain Kohn- 
Sham like equations by the variational method, 

dE LDA + G {p) _^ dE LDA +G{p) ^ 

d^a dmr 
The charge density under GA can be constructed as: 
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In practice, two steps are followed. (1) For fixed mr, 
the |*o) is optimized by solving the effective Kohn-Sham 
equations. This step is basically the same as all other 
LDA calculations, with additional little cost for the pro- 
jection to local orbitals; (2) For fixed l^o), mr are opti- 
mized by solving set of linear equations 11( , which have 
dimension of number of configurations T (for instance 




FIG. 1: The calculated band structure and density of states 
for SrVC>3 using LDA and LDA+Gutzwiller method. 
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FIG. 2: The calculated total energy as function of volume for 
different phases of bulk Fe by LDA+Gutzwiller method. The 
bcc FM ground state is correctly predicted. 



2 10 for five d orbitals, and the number of configurations 
can be reduced significantly by considering the symme- 
try). The second step is addition to the LDA or LDA+U 
scheme, however, the computational cost is rather small 
compared with the cost by wave function optimization 
(i.e., the first step), because only set of linear equations 
needed to be solved [llj ■ The final solution of the system 
is obtained until the charge density (and total energy) is 
self-consistently converged. 

The above proposed LDA+Gutzwiller method was im- 
plemented in our BSTATE (Beijing Simulational Tool 
for Atom Technology) code [13j], which uses plane- wave 
ultra-soft pseudo-potential method. For the local basis, 
the physical choice is the wannier function constructed 
from the atomic orbital. This was done by using the pro- 
jected wannier function method [l2T | in our calculations. 
For the interaction term, the on-site density-density in- 
teractions are considered in the present studies. Several 
typical systems have been calculated, and the results are 
discussed below. 

1. Non-magnetic correlated metal: SrVOz 

SrVC>3 is a correlated metal with Sd-t^g configura- 
tion. The simple cubic perovskite crystal structure and 
non-magnetic electronic state make it an idea test mate- 
rial 14]. The one-particle spectrum calculated with LDA 
is poorly compared with experiments [l5| . The main 
problems are: (1) the calculated band width is about 40% 
wider than photoemission observation; (2) the estimated 
effective mass is about 2-3 times lower than experimental 
results from specific heat and susceptibility; (3) the pho- 
toemission peak observed around 2eV below Fermi level 
(Ep) is not understood. 

However, most of the features can be understood from 
our LDA+Gutzwiller calculations. We choose the 3d 
Wannier function as the interacting local orbits with ef- 
fective interaction energy £7=5. OeV following the liter- 
atures Fig.l shows the calculated band structure 
and the density of states of quasi-particle spectrum. It 
is clear that the band width is reduced by about 40% 
compared with LDA. From the calculated quasi-particle 



spectrum it is straight forward to calculate the effective 
mass enhancement (rn*/m), which is 2.1 times larger 
than LDA. All the results are in good consistency with 
experiments [TBI ], and much better than that obtained 
from LDA. Finally, it is worth to make the comment: 
since all the quasi particle part has been well treated in 
the present scheme (incoherent part is not included yet), 
and no DOS is found around -2eV region, it further sug- 
gests that the observed photoemission peak in this region 
should be incoherent as suggested by other studies [14 1 . 
2. Magnetic correlated metal: Fe and Ni 

Bulk Fe and Ni are typical magnetic metal with inter- 
mediate correlations. For Fe, the LDA fails to predict 
the bcc FM ground state, although GGA correctly do 
so. Both LDA and GGA overestimate the band width 
by about 10-20% compared with experiements [l6j]. For 
the Ni, the problems are more serious, the band width is 
overestimated by about 30%, and the spin polarization 
is hardly compared with experiments [9|. For such in- 
termediately correlated metal, the LDA+U method may 
improve one or two discrepancies to some extent, while 
the price to be payed is that other predicted properties 
become even worse than LDA. 

For our calculations, we first determine the effective U 
from constrained-LDA calculations to be 7.0eV (Fe) and 
9.0eV (Ni) for the choice of our local orbital, and Hund's 
coupling J is fixed to be l.OeV (a common choice for 3d 
elemental metal). The results, summarized in Table I 
and Fig. 2, show that most of the discrepancies are sys- 
tematically corrected compared with experiments, sug- 
gesting the advantages of present scheme. In particular, 
the following improvements are significant: (1) the band 
width rcnormalization is correctly predicted and effective 
mass enhancement (specific heat coefficient) can be well 
compared with experiments. (2) the calculated lattice 
parameter, bulk modulus and magnetic moment show 
systematic improvement. Finally, we want to point out 
that for the down-spin band of Ni, around the X point 
near Ep, one of the bands predicted by LDA to be above 
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FIG. 3: The calculated band structure of fee FM Ni. 



TABLE I: The calculated property parameters for bec FM 
Fe and fee FM Ni in comparison with experimental results. 
They are equilibrium lattice constant aa, bulk modulus B, 
spin magnetic moment M, specific heat coefficient 7, and the 
occupied energy band width W. The experimental data are 
from Ref. [13] ■ 
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W(eV) 




LDA 


5.21 


227 


2.08 


2.25 


3.6 


Fe 


LDA+G 


5.39 


160 


2.30 


3.52 


3.2 




Exp. 


5.42 


168 


2.22 


3.1,3.69 


3.3 




LDA 


6.49 


250 


0.59 


4.53 


4.5 


Ni 


LDA+G 


6.61 


188 


0.50 


6.9 


3.2 




Exp. 


6.65 


186 


0.42,0.61 


7.02 


3.2 



Ef, is now correctly predicted to be below Ef, in con- 
sistency with recent AREPS data 
3. AF correlated Insulator: NiO 

In the present LDA+Gutzwiller scheme, both effec- 
tive on-site level renormalization and kinetic renormal- 
ization are included, while only the former is considered 
in LDA+U scheme. For the intermediately correlated 
metallic systems, the kinetic energy renormalization is 
significant, this is reason why LDA+U scheme fails. How- 
ever, for the large U limit, where effective on-site level 
renormalization dominate, such as the AF correlated in- 
sulator NiO, our calculations suggest that the kinetic 
energy renormalization factor q a is almost utility, and 
the obtained electronic structure from LDA+Gutzwiller 
method is almost the same as that of LDA+U. How- 
ever, it should be noted that this is true only for the 
cases where q a close to 1. For some of the AF insula- 
tors, the localized moment is far away from utility, we 
expect that the kinetic energy renormalization will also 
contribute. In these cases, the present LDA+Gutzwiller 
approach will give a better results than LDA+U even for 
the AF ordered insulators. We will leave this issue for 
future studies. 

In summary, we have shown that the Gutzwiller ap- 
proach can be well combined with the DFT through vari- 
ational principal. As the results, a fully self-consistent 
LDA+Gutzwiller method is developed. The method not 
only keeps the accuracy and parameter-free character of 
LDA-type calculations (such as the total energy calcu- 
lations), but also automatically cover the region from 



weakly to strongly correlated systems. The quasi particle 
spectrum is properly described by taken into account the 
kinetic energy renormalization, and the calculated results 
for several typical systems demonstrate the strong advan- 
tage over LDA and LDA+U results. On the other hand, 
the scheme itself is easy to be implemented and compu- 
tationally much cheaper than the LDA+DMFT method. 
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